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Abstract — This paper presents novel adaptive space-time 
reduced-rank interference suppression least squares algorithms 
based on joint iterative optimization of parameter vectors. The 
proposed space-time reduced-rank scheme consists of a joint 
iterative optimization of a projection matrix that performs 
dimensionality reduction and an adaptive reduced-rank 
parameter vector that yields the symbol estimates. The proposed 
techniques do not require singular value decomposition (SVD) 
and automatically find the best set of basis for reduced-rank 
processing. We present least squares (LS) expressions for the 
design of the projection matrix and the reduced-rank parameter 
vector and we conduct an analysis of the convergence properties 
of the LS algorithms. We then develop recursive least squares 
(RLS) adaptive algorithms for their computationally efficient 
estimation and an algorithm for automatically adjusting the 
rank of the proposed scheme. A convexity analysis of the LS 
algorithms is carried out along with the development of a 
proof of convergence for the proposed algorithms. Simulations 
for a space-time interference suppression application with a 
DS-CDMA system show that the proposed scheme outperforms 
in convergence and tracking the state-of-the-art reduced-rank 
schemes at a comparable complexity. 

Index Terms — pace-time adaptive processing, interference sup- 
pression, spread spectrum systems, iterative methods, least 
squares algorithms, pace-time adaptive processing, interference 
suppression, spread spectrum systems, iterative methods, least 
squares algorithms, s 



I. Introduction 

SPACE-TIME adaptive processing (STAP) techniques have 
become a fundamental enabling technology of modern 
systems encountered in communications |fl~), radar and sonar 
El, 0, and navigation [4]. The basic idea is to gather data 
samples from an antenna array and process them both spatially 
and temporally via a linear combination of adaptive weights. In 
particular, STAP algorithms have found numerous applications 
in modern wireless communications based on spread spectrum 
systems and code-division multiple access (CDMA) 0, |6|. 
These systems implemented with direct sequence (DS) sig- 
nalling are found in third-generation cellular telephony |f71, |8l, 
|9l , indoor wireless networks iflOl , satellite communications, 
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ultra-wideband technology ifTTl and are being considered for 
future systems with multi-carrier versions such as MC-CDMA 
and MC-DS-CDMA Q2|, and in conjunction with multiple 
antennas (TJ). The advantages of spread spectrum systems 
include good performance in multi-path channels, flexibility 
in the allocation of channels, increased capacity in bursty 
and fading environments and the ability to share bandwidth 
with narrowband communication systems without performance 
degradation [5|. 

There are numerous algorithms with different trade-offs be- 
tween performance and complexity for designing STAP tech- 
niques 1(141 . Among them, least squares (LS)-based algorithms 
are often the preferred choice with respect to convergence 
performance. However, when the number of filter elements 
in the STAP algorithm is large they require a large number of 
samples to reach its steady-state behavior and may encounter 
problems in tracking the desired signal. Reduced-rank STAP 
techniques lfl5l - ll4"0l are powerful and effective approaches 
in low-sample support situations and in problems with large 
filters. These algorithms can effectively exploit the low-rank 
nature of signals that are found in spread spectrum com- 
munications. Their advantages are faster convergence speed 
and better tracking performance than full-rank techniques 
when dealing with a large number of weights. It is well 
known that the optimal reduced-rank approach is based on the 
singular value decomposition (SVD) of the known input data 
co variance matrix R [15). However, this co variance matrix 
must be estimated. The approach taken to estimate R and 
perform dimensionality reduction is of central importance and 
plays a key role in the performance of the system. Numerous 
reduced-rank strategies have been proposed in the last two 
decades. Among the first methods are those based on the 
SVD of time-averaged estimates of R lfT5l - l20l . in which 
the dimensionality reduction is carried out by a projection 
matrix formed by appropriately selected eigenvectors com- 
puted with the SVD. An effective approach to address the 
problem of selection of eigenvectors, known as the cross- 
spectral method, and that results in improved performance was 
considered in lETl . Iterative algorithms that avoid the SVD but 
do not fully exploit the structure of the data for reduced-rank 
processing were reported in El . Il23l . A more recent and 
elegant approach to the problem was taken with the advent 
of the multistage Wiener filter (MSWF) J24], which was later 
extended to adaptive versions by Honig and Goldstein in l25l , 
STAP applications l26l and other related techniques ll27l . 
Another method that was reported about the same time as 
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the MSWF is the auxiliary vector filtering (AVF) algorithm 
l|28l - ll32l . A reduced-rank method based on interpolated filters 
with time- varying interpolators was reported in ||34l , ll35l , IT361 
for temporal processing and an associated STAP version was 
considered in [37), however, this approach shows significant 
performance degradation with small ranks. A key limitation 
with the existing reduced-rank STAP techniques is the lack 
or a deficiency with the exchange of information between the 
projection matrix that carries out dimensionality reduction and 
the subsequent reduced-rank filtering. 

In this work we propose reduced-rank STAP LS algorithms 
for interference suppression in spread spectrum systems. The 
proposed algorithms do not require SVD and prior knowledge 
of the reduced model order. The proposed reduced-rank STAP 
scheme consists of a joint iterative optimization of a projection 
matrix that performs dimensionality reduction and is followed 
by an adaptive reduced-rank filter. The key aspect of the 
proposed approach is to exchange information between the 
tasks of dimensionality reduction and reduced-rank processing. 
The proposed STAP scheme builds on the temporal scheme 
first reported in 11381 with stochastic gradient algorithms and 
extends it to the case of spatio-temporal processing and to 
a deterministic exponentially-weighted least squares design 
criterion. We develop least squares (LS) optimization algo- 
rithms and expressions for the joint design of the projection 
matrix and the reduced-rank filter. We derive recursive LS 
(RLS) adaptive algorithms for their computationally efficient 
implementation along with a complexity study of the proposed 
and existing algorithms. We also devise an algorithm for 
automatically adjusting the rank of the filters utilized in the 
proposed STAP scheme. A convexity analysis of the proposed 
LS optimization of the filters is conducted, and an analysis of 
the convergence of the proposed RLS algorithms is also carried 
out. The performance of the proposed scheme is assessed 
via simulations for a space-time interference suppression ap- 
plication in DS-CDMA systems. The main contributions of 
this work are summarized as follows: 1) A reduced-rank 
STAP scheme for spatio-temporal processing of signals; 2) 
LS expressions and recursive algorithms for STAP parameter 
estimation; 3) An algorithm for automatically adjusting the 
rank of the filters; 4) Convexity analysis and convergence proof 
of the proposed LS-based algorithms. 

This work is organized as follows. Section II presents 
the space-time system model, and Section III states the 
reduced-rank estimation problem. Section IV presents the 
novel reduced-rank scheme, the joint iterative optimization and 
the LS design of the filters. Section V derives the RLS and 
the rank adaptation algorithms for implementing the proposed 
scheme. Section VI develops the analysis of the proposed 
algorithms. Section VII shows and discusses the simulations, 
while Section VIII gives the conclusions. 

II. Space-Time System Model 

We consider the uplink of DS-CDMA system with symbol 
interval T, chip period T c , spreading gain N = T/Tc, K 
users, multipath channels with L propagation paths and 
L < N. The system is equipped with an antenna that consists 



of a uniform linear array (ULA) and J sensor elements Q, 0. 
In the model adopted, the intersymbol interference (ISI) span 
and contribution are functions of the processing gain N and L 
@. For instance, we assume that L < N which results in the 
interference between 3 symbols in total, the current one, the 
previous and the successive symbols. The spacing between the 
ULA elements is d = A c /2, where A c is carrier wavelength. 
We assume that the channel is constant during each symbol, 
the base station receiver is perfectly synchronized and the 
delays of the propagation paths are multiples of the chip rate. 
The received signal after filtering by a chip-pulse matched 
filter and sampled at the chip period yields the JM x 1 received 
vector at time i 

K 

r\i] =Y, A kh[i - l]p k [i - 1] + A k b k \i\p k [i] 

fc=i *• ' 

+ A k b k [i + l]p k [i + 1} + n[i], 

where M = N + L — 1, the complex Gaussian noise vector 
is n[i] = [ni[i] ... njM{i]] T with i?[n[i]n H [i]] = a 2 I, 
(•) T and (•) denote transpose and Hermitian transpose, 
respectively, and E[-] stands for expected value. The spatial 
signatures for previous, current and future data symbols are 

p k [i-i} = F k n k [i-i}, 

p k \i}=F k H k \i], (2) 
p k [i + 1} = T k n k [i + 1], 

where f k , T k and T k are block diagonal matrices 
with versions of segments of the signature sequence s k = 
[cifc(l) . . . a k (N)] T of user k shifted down by one position 
(one-chip) and given by T k = diag(c fc , C k , . . . , C k j , T k = 

diag(Cfc,C fc ,...,C fc ), and T k = diag(c fe , C k , . . . , C fe ) . 

The structure of the M xL matrices C k , C k and C k is detailed 
in 0. The JL x 1 space-time channel vector is given by 

n k [i] = [hlM hlM ... |h£ J _ 1 [i]] T , (3) 

where h k j[i] = [/i^oIX] ■ • ■ h k l \_ l \i]] T is the L x 1 channel 
vector of user k at antenna element I with their associated 
directions of arrival (DoAs) <fi k , m . The DoAs are assumed 
different for each user and path 11241 . 

III. Reduced-Rank STAP for Interference 
Suppression and Problem Statement 

In this section, we outline the main problem of STAP 
design for interference suppression in spread spectrum systems 
and we consider the design of reduced-rank STAP algorithms 
using a least squares approach. The main goal of the STAP 
algorithms is to jointly perform temporal filtering with spatial 
filtering (beamforming) through adaptive combination of filter 
coefficients. 

Let us consider the space-time received signals of the 
previous section and the data organized in JM x 1 vectors 
r[i]. In order to process this data vector, one can design a 
STAP algorithm that consists of a JM x 1 filter w[i] = 
[w® «4 • ■ • w ja/] T ' which adaptively and linearly combines 
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its coefficients with the received data samples to yield an esti- 
mate x[i] = w H [i]r[i]. The design of w[i] can be performed 
via the minimization of the exponentially weighted LS cost 
function 



C(w[i]) = ^2\ i - l \d[l\-w B [%\r\ 



(4) 



1=1 



where d[l] is the desired signal and A stands for the forgetting 
factor. Solving for w[i], we obtain 



w\i\ 



fl _1 [*>[i] 



(5) 



where R[i] = Yl]=i ^ r[l]r H [I] is the time-averaged cor- 
relation matrix and p[i] = 22; =1 A d*[7]r[Z]] is the cross- 
correlation vector. 

A problem with STAP algorithms is that the laws that 
govern their convergence and tracking behavior imply that the 
performance is a function of JM, the number of elements 
in the filter. Thus, large JM implies slow convergence and 
poor tracking performance. A reduced-rank STAP algorithm 
attempts to circumvent this limitation by exploiting the low- 
rank nature of spread spectrum systems and performing dimen- 
sionality reduction. This dimensionality reduction reduces the 
number of adaptive coefficients and extracts the key features 
of the processed data. It is accomplished by projecting the re- 
ceived vectors onto a lower dimensional subspace. Specifically, 
consider a JM x D projection matrix To [i] which carries out 
a dimensionality reduction on the received data as given by 



(6) 



where in what follows all Z?-dimensional quantities are de- 
noted with a "bar." The resulting projected received vector 
f\i\ is the input to a tapped-delay line represented by f " 



vector w[i\ = [w 
output is 



.r ?. 



The reduced-rank 



w H \i]r\i] 



If we consider the LS design in with the reduced 
parameters we obtain 



R l [i]p[i] 



^[i\R\i}Tn[i\ 



where R[i] = ELi A f M = T i 
reduced-rank correlation matrix, p[i] = YTi=i A l_! d*[Z]i 
T^[i]p[i] is the cross-correlation vector of the reduced-rank 
model. The associated sum of error squares (SES) for a rank- 
le STAP is expressed by 



SES = o d — p [i]R \i]p[i] 

= o 2 d -p H \ l ]T D ^{T H D [m\i]TD[i])- 1 T H D m^ 



(9) 



Ez=i ^ 'M(0| 2 - The development above shows 



where a d 

us that the key aspect for constructing reduced-rank STAP 
schemes is the design of To [i] since the SES in (0 depends 
on p[i], R[i] and Tu[ij. The quantities p[i] and R[i] are 
common to both reduced-rank and full-rank STAP designs, 
however, the projection matrix To[i] plays a key role in the 
dimensionality reduction and in the performance. The strategy 
is to find the most appropriate trade-off between the model 
bias and variance ifTSl by adjusting the rank D and exchanging 



information between To[i] and w[i]. For instance, evaluating 
numerically the SES expression in equation (0 one can verify 
the convergence and steady state performance of reduced-rank 
STAP algorithms. Next, we present the proposed reduced-rank 
STAP approach. 

IV. Proposed Reduced-Rank STAP and Least 
Squares Design 

In this section we detail the principles of the proposed 
reduced-rank STAP scheme and present a least squares (LS) 
design approach for the filters. The proposed reduced-rank 
STAP scheme is depicted in Fig. 1 and is formed by a 
projection matrix Td[z] with dimensions JM x D that is 
responsible for the dimensionality reduction and a D X 1 
reduced-rank filter w[i]. The JM x 1 received data vector r[i] 
is mapped by T o[i] into a D X 1 reduced-rank data vector 
r[i). The reduced-rank filter w[i] linearly combines r[i) in 
order to yield a scalar estimate x[i\. The key strategy of the 
proposed framework lies in the joint design of the projection 
matrix To[i] and the reduced-rank filter w[i] according to the 
LS criterion. The exchange of information between Tjj[i] 
and w[i] is different from the MSWF J24]-||26) and the AVF 
techniques ||28l - ||32l . In particular, the expressions of the 
filters obtained for the proposed reduced-rank STAP scheme 
allow a more efficient introduction of the bias than that of the 
MSWF and the AVF by alternating the recursions for Td[i\ 
and w[i]. In addition, the proposed STAP scheme is based on a 
subspace projection designed according to a joint and iterative 
minimization of the LS cost function and which achieves better 
performance than the Krylov subspace of the MSWF and the 
AVF. 



Projection Matrix 



Reduced-Rank Filter 
m [ij 



l[iJ=W (i|r(i] 



JM xl 
I 



' (') 



Design Algorithm -4— 



Fig. 1. Proposed Reduced-Rank STAP Scheme. 

Let us now detail the quantities involved in the proposed 
reduced-rank STAP scheme. Specifically, the projection matrix 
To[i] is structured as a bank of D full-rank filters with 
dimensions JM x 1 which are described by 

iT 



L JM.d\ ' " 



1, . .., D. (10) 
The filters td[i] are then gathered and organized, yielding 



T D [i\=[tf\tf 



t [l] 1 



(ID 



The output estimate x[i] of the reduced-rank STAP scheme 
can be expressed as a function of the received data r[i], the 
projection matrix To[i] and the reduced-rank filter w[i] as 
given by 



w H \i}T^\i]r\i] =w H \i]f\i}. 



(12) 
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Interestingly, for D = 1, the proposed STAP scheme becomes 
a conventional full-rank STAP algorithm with an additional 
weight parameter wd that can be seen as a gain. For D > 1, 
the signal processing tasks are changed and the full-rank filters 
td[i] perform dimensionality reduction and the reduced-rank 
filter estimates the desired signal. 

In order to design the projection matrix Tu[i] and the 
reduced-rank filter w[i] we need to adopt an appropriate design 
criterion. We will resort to an exponentially-weighted LS ap- 
proach since it is mathematically tractable and results in joint 
optimization algorithms that can track time-varying signals by 
adjusting the forgetting factor A. The design of the proposed 
scheme amounts to solving the following optimization problem 

i 

[T Aopt [i],™? pt [i]] = arg min £ A 4 "'^] - w H [{\T^ 

(13) 

In order to solve the above minimization problem, let us then 
consider the cost function 

i 

c(T D ii],w H w) = Y,^ l \dii]-w H mzmi}\ 2 ], (H) 

i=i 

Minimizing ( TT4l with respect to [i] we obtain 

T D , opt \i] = R-^Pd^R^Ii}, (15) 

where Pfl[i] = ^i=i X l ~ l d*[l]r[l]w H [i], the time-averaged 
correlation matrix is R[i] = 2;=i X l ~ l r[l]r H [I] and Rw[i] = 
w[i]w H [i]. Note that we have opted for computing 
as Jie,[i] = J2]=i X w[l]w [I] with a regularization term 
introduced at the beginning of the iterations in order to 
allow the computation of its inverse. For this reason, the 
latter approach will be adopted for the derivation of adaptive 
algorithms. Minimizing ( TBI with respect to w[i], the reduced- 
rank filter becomes 

w opt [i] = -R -1 [i]p[i], (16) 
where p\i] rg )0pt [i] E|=i A*-'d*[i]r[/] 

Tl, opt V\Y!i=i^ l r[l]r H [l]T D ^\i]. The associated 
SES for the proposed reduced-rank STAP scheme is 



SES = ffj-p w [i]JR \i]p[i\ 



(17) 



where a 2 d = Ya=i ^ \d[l]\ 2 ■ Note that the expressions in 
(fT3T > and (\M are not closed-form solutions for tu ptH an d 
Td, opt [i] since ( TT~5b is a function of «jopt[*] an d *0 depends 
on and To op t[i]. Therefore they have to be iterated with an 
initial guess to obtain a solution. The expressions in (fTBT l and 
< fT~6b require the inversion of matrices, which entails cubic 
complexity with JM and D. Computing the SES in ( TP7I ). it 
can be numerically verified the convergence and steady state 
performances of reduced-rank STAP algorithms, namely, the 
proposed, the MSWF |25] and the AVF J32). In order to reduce 
the complexity, we will derive RLS algorithms in the next 
section. The rank D must be set by the designer to ensure 
appropriate performance or a mechanism for automatically 
adjusting the rank should be adopted. We will also present 
an automatic rank selection algorithm in what follows. 



V. Proposed RLS and Rank Selection Algorithms 

In this section we propose RLS algorithms for efficiently 
implementing the LS design of the previous section and 
estimating the filters Tr}, pt[i) an d ^opt[i] with the filters 
Td[g'] and w[i], respectively. We also develop rank selection 
algorithms for automatically adjusting the rank D of the 
proposed STAP algorithm. An analysis of the computational 
requirements of the proposed and analyzed algorithms is also 
included. 



A. Proposed RLS Algorithm 

In order to derive an RLS algorithm for the proposed 
scheme, we consider ( fTst and derive a recursive procedure 
iJb^^Sjijputing the parameters of To[i]. Let us define 

p[i] = R- 1 [i\, 

Q^] = R^[i-l], (18) 
P D [%\ = XP D [i - 1] + d* [»] r [i]w H [i] , 

and rewrite the expression in dl5t as follows 

T D [i]=P[ l ]P D [ l ]Q a [i] 

= XP\i]P D [i - l]Q^,[i]+d*[{\P\i]r\i]w H [i]Q^\i] 
= T D [i-l}- k[i\P[i - l]P D [i - l]Qa,]i] 

+ d*[i]P[i]r[i]w H [i]Q & [i] 
= T D [i-l}- k[i]P[i - l]P D [i - l]Q a \i] 

+ d*\i]k\i]w H [i]Q^[i] 

(19) 

By defining the vector t[i] = Q ffi [i]ti)[i] and using the fact 
that f H [i - 1] = r H [i - l]T D [i - 1] we arrive at 

T D [i\ =T D [i-l] + k[i\(d*\i}t H [i\ - r H [i]), (20) 

where the Kalman gain vector for the computation of T o[i] 
is 

A- x P[i - l\r[i] 



k\i] 



1 + \- x r H \i}P[i - l]r[*] 



(21) 



and the update for the matrix P[i] employs the matrix inver- 
sion lemma lfT4l 

P[i] = \- l P[i - 1] - \- l k\i}r H \i}P[i - 1] (22) 

the vector t[i] is updated as follows 

X^Qvli-lWi-l] 



(23) 



1 + A- 1 ™^ - l]Q a [i - l]w[i - 1] 

and the matrix inversion lemma is used to update Q ffl [i] as 
described by 

Qw\i] = X^Q^i - 1] - \-H[i]w H [i - l]g e [i - 1], (24) 

The equations (f20b - (124-b constitute the first part of the pro- 
posed RLS algorithm and are responsible for calculating the 
projection matrix Td[z]. 

In order to derive a recursive update equation for the 
reduced-rank filter w[i], we consider the expression in 
< n~6] > with its associated quantities, i.e., the matrix R[i] = 



5 



(25) 



ELi y~ l r[l]f H [l] and the vector p[i) = £j =1 ^~ l d*[l]r[ 
Let us define 

*[i]=fl- 1 [i], 
p[i] = Xp[i - 1] + d*[i]f [i], 
and then we can rewrite ( [TBI in the following alternative form 

w[i] = w[i - l]+k[i][d*[i] -f H [i]w[i - 1]] (26) 

By defining = d[i] — w H [i — l]f H [i] we arrive at the 
proposed RLS algorithm for obtaining w[i] 

w[i] = w[i-l] + k[i]C[i\: (27) 
where the so-called Kalman gain vector is given by 

A _1 #[« - l]f[i] 



ki] 



1 + A- 1 ^^']*^ - l]f[i] 



(28) 



and the update for the matrix inverse 3? [i] employs the matrix 
inversion lemma |[T4l 

= A -1 *[i - 1] - A _1 fc[i]f ff [i]$[i - 1]. (29) 

It should be noted that the proposed RLS algorithm given in 
(|27]>-(|291 is similar to the conventional RLS algorithm IfBl . 
except that it works in a reduced-rank model with a D x 1 
input r[i] = Tp[i]r[i], where the JM x D matrix Td is the 
projection matrix responsible for dimensionality reduction. 

B. Rank Selection Algorithm 

The performance of the RLS algorithm described in the 
previous subsection depends on the rank D. This motivates the 
development of methods to automatically adjust D on the basis 
of the cost function. Unlike prior methods for rank selection 
which utilize MSWF-based algorithms 11251 or the cross- 
validation approach used with AVF-based recursions 11321 . we 
focus on an approach that determines D based on the LS 
criterion computed by the filters Tjj[i] and w^ D > [i], where the 
superscript ( D > denotes the rank used for the adaptation. Al- 
though there are similarities between the algorithm described 
here and the one reported in ||2"5l , the algorithm presented 
here differs from ll2~5l in that it clearly details the strategy for 
updating Td[z] and w^ D \i], defines the maximum (-D ma x) 
and minimum (D mm ) values for the rank D allowed and works 
with extended filters for reduced complexity. The method for 
automatically selecting the rank of the algorithm is based on 
the exponentially weighted a posteriori least-squares type cost 
function described by 

i 

C ap (T D [z],*( D )[z]) = ^a^|#]-^ ^[i]T D [i]r[l]\^ 

i=i 

(30) 

where a is the forgetting factor and w^[i] is the reduced- 
rank filter with rank D. For each time interval i, we can select 
D which minimizes C(Tjj[i], w^ D '[i]) and the exponential 
weighting factor a is required as the optimal rank varies as 
a function of the data record. The dimensions of Td[«] and 
w^li] are extended to M x £) max and -D max , respectively, 
and the associated matrices R[i], Pj[i] and Q^[i] should 
be compatible for adaptation. Our strategy is to consider the 



adaptation with the maximum allowed rank Z? max and then 
perform a search with the aim of finding the best rank within 
the range D min to £> max . To this end, we define Tu[i] and 



,(£>) 



as follows: 



T D [i] = [ ti[i] ... ti min [i] 



(31) 



w^li] = [ wi[i] ... w Dmia [i] ... wn mal [i] J 
The proposed rank selection algorithm is given by 

A> pt [i] = arg min C ap (T D [i\,w^[i\), (32) 

where d is an integer, D m i n and D max are the minimum and 
maximum ranks allowed for the reduced-rank filter, respec- 
tively. Note that a smaller rank may provide faster adaptation 
during the initial stages of the estimation procedure and a 
greater rank usually yields a better steady-state performance. 
Our studies reveal that the range for which the rank D 
of the proposed algorithms have a positive impact on the 
performance of the algorithms is limited, being from D m i n = 3 
to -D max = 8 for the reduced-rank filter recursions. These 
values are rather insensitive to the system load (number of 
users), to the number of array elements and work very well 
for all scenarios and algorithms examined. The computational 
complexity of the proposed rank selection algorithm with 
extended filters is equivalent to the computation of the cost 
function in ( f30b and requires 3(Z? max — Z? m ; n ) + 1 additions and 
a sorting algorithm to find the best rank according to (l32l . An 
alternative strategy to using extended filters is the deployment 
of multiple filters with the rank selection algorithm in (l32t 
that determines the best set of filters for each time interval. 
Specifically, this approach employs Z3 max — D m [ n + 1 pairs of 
filters and has a very high complexity. 

A second approach that can be used is a mechanism based 
on the observation of the columns of Tjj[i] and a stopping 
rule, as reported in ||251 . The method performs the following 
optimization 



D, 



opt [ 



arg max 

d 



\\td[i\\ 



>s, 



(33) 



where Px d (x) is the orthogonal projection of the vector x onto 
the subspace Td and 5 is a small positive constant. In l25ll . it 
has not been discussed the use of a range of values for allowing 
the selection, however, we found that it is beneficial in terms 
of complexity to restrict the optimization to an appropriate 
range of values Z? max to D m \ n as with the previous method. 

Another possibility for rank selection is the use of the cross- 
validation (CV) method reported in 11321 . This approach selects 
the filters' lengths which minimize a cost function that is 
estimated based on observations (training data) that have not 
been used in the process of building the filters themselves as 
described by 



C C v(T D W , w (D) [»])=£ a*" 1 \d{l)-w^ D) [i]T A(i/I) [i]r[ 

(34) 



We consider here the same "leave one out" approach as in 
For a given data record of size i, the CV approach chooses 
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TABLE I 

Computational complexity of RLS algorithms. 



TABLE II 

Computational complexity of the proposed rank 
selection algorithm with multiple filters. 



Algorithm 
Full-rank Q4) 

Proposed 
MSWF (25) 
AVF m 



Additions 

3(JM) 2 - 2JM + 3 

3(JM) 2 - 2.JM + 3 
6D 2 -8D + 3 

D{JM) 2 + (JM) 2 + 6D 2 
-8D + 2 

D({JM) 2 +3(JM- l) 2 ) 
+D(5(JM - 1) + 1) 
2JM - 1 



Multiplications 

6( JM) 2 + 2JM + 2 

Proposed with 
7( JM) 2 + 2JM Multiple Filters 

7D 2 + 9D 



2(Anax - Anin) + 1 

(D max -Dfain + l)x 

(3( JM) 2 - 2JM + 3 
+6ALx - 8D max + 3) 



(Anax ~ Anin + 1)X 

(7( JM) 2 + 2JM 



D(JM) 2 + (JM) 2 TABLE III 

2DJM + 3D + ^Computational complexity of remaining rank selection 

ALGORITHMS. 

D(4(JM) 2 + 4JM + 1) 
4JAf + 2 



the filter io^' [i] that performs the following optimization 

£> opt [i] = axg min CcvC^], tj)M[i]), (35) 
de{i,2,...} 

The main difference between this and the other algorithms 
presented lies in the use of CV, which leaves one sample out in 
the process, and the use of the constraint on the allowed filter 
lengths. In the simulations, we will compare the rank selection 
algorithms and discuss their advantages and disadvantages. 

C. Computational Complexity 

In this part of the work, we detail the computational 
complexity requirements of the proposed RLS algorithms 
and compare them with those of existing algorithms. We 
also provide the computation complexity of the proposed 
and existing rank selection algorithms. The computational 
complexity expressed in terms of additions and multiplications 
is depicted in Table I for the RLS algorithms, the complexity 
of the proposed rank selection algorithm with multiple filters 
including the proposed RLS algorithm is illustrated in Table 
II, and that of the remaining rank selection techniques is given 
in Table III. 




O Proposed 

• • MSWF 

• AVF 
Full-Rank 



20 40 60 80 100 120 
JM 



20 40 60 80 100 120 
JM 



Fig. 2. Complexity in terms of additions and multiplications against 
number of input samples (JM) and D — 4 . 



In the case of the proposed reduced-rank RLS algorithm 
the complexity is quadratic with (JM) 2 and D 2 . This corre- 



Algorithm 
Prop osed with 
Extended Filter 



Projection with 
Stopping Rule (25) 

CV El 



Additions 

2(D 

max 

-D 

mm j 



+ 1 



2(2 JM - l)x 
(D max -D mm 

(2JM - l)x 

(2(D 

max D m i n 

) + i) 



Multiplications 



((JM) 2 + JM + 1) x 

max 

-D m i n 

+ i) 



+ l)x 

JM + 1 



sponds to a complexity slightly higher than the one observed 
for the full-rank RLS algorithm, provided D is significantly 
smaller than JM, and comparable to the cost of the MSWF- 
RLS J25] and the AVF [32). In order to illustrate the main 
trends in what concerns the complexity of the proposed and 
analyzed algorithms, we show in Fig. |2]the complexity against 
the number of input samples JM. The curves indicate that 
the proposed reduced-rank RLS algorithm has a complexity 
lower than the MSWF-RLS algorithm |25j and the AVF l32l . 
whereas it remains at the same level of the full-rank RLS 
algorithm. 

The proposed rank selection algorithm with multiple filters 
has a number of arithmetic operations that is substantially 
higher than the other compared methods since it requires the 
computation of D max — D m - in + 1 pairs of filters with the 
proposed RLS algorithms simultaneously. We show the overall 
cost of this algorithm separately in Table II. The computational 
complexity of the remaining rank selection algorithms includ- 
ing the proposed and the existing rank selection algorithms 
is shown in Table III. From Table III, we can notice that 
the proposed rank selection algorithm with extended filters 
is significantly less complex than the existing methods based 
on projection with stopping rule [|25j and the CV approach 
ll32l . Specifically, the proposed rank selection algorithm with 
extended filters only requires 2(D max — D m - m ) additions, as 
depicted in the first row of Table III. To this cost we must add 
the operations required by the proposed RLS algorithm, whose 
complexity is shown in the second row of Table I using D max 
according to the procedure outlined in the previous subsection. 
The complexities of the MSWF and the AVF algorithms are 
detailed in the third and fourth rows of Table I. For their 
operation with rank selection algorithms, a designer must add 
their complexities in Table I to the complexity of the rank 
selection algorithm of interest, as shown in Table III. 

VI. Analysis 

In this section, we conduct a convexity analysis of the 
proposed optimization that is responsible for designing the 
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filters To[i] and w[i] of the proposed scheme. We show that 
the proposed optimization leads to a problem with multiple 
solutions, and we discuss the properties of the method. In par- 
ticular, we conjecture that it leads to a problem with multiple 
and possibly identical minimum points. This is corroborated 
by numerous studies that verify that the method is insensitive 
to different initializations (except for the case when Tfl[i] is 
a null matrix and which annihilates the received signal) and 
that is always converge to the same point of minimum. We 
also establish the convergence of the proposed optimization 
algorithm, showing that the sequence of filters T,o[i] and w[i] 
produces a sequence of outputs that is bounded. 

A. Convexity Analysis of the Proposed Method 

In this part, we carry out a convexity analysis of the 
proposed reduced-rank scheme and LS optimization algorithm. 
Our approach is based on expressing the output of the pro- 
posed scheme in a convenient form that renders itself to 
analysis. Let us consider the proposed optimization method 
in ( [T3b and express it by an expanded cost function 



where U[i] is a D(JM + 1) x D(JM + 1) matrix which 
contains !R[i]. At this stage, we can alternatively express the 
cost function in (f36t as 



(39) 



C{q\i])=Y,\d[l]-q H mm 



i=i 



We can examine the convexity of the above by computing the 
Hessian (i?)with respect to q[i] using the expression 

d d(C(q[i}) 



H 



(40) 



c(T D \i],w H \i}) =j2*- l \d[i]-™ H m%mi}\ 2 ] 
i=i 



1=1 



1=1 



Y^^- l d[l]r H [l]T D [i}w[i} 
i=i 



E 

i=i 



\ i ~ l w H [i]rf [i]r[l]r H [l]T D [i]w[i] , 

(36) 



w H [i\ 


' r[i] . 
r[i] 


.. 

■. 


T 


- 8 $[i\ - 




... . 


. . r[i] _ 







In order to proceed, let us express x[i] in an alternative and 
more convenient form as 

D 

x[i\ = w H \i]T H D \i]v[^ = w H [i\^2T^\i}r[i]v d 



(37) 



= w H [i\M T [i\ S * v \i] 

where is a DJM x D block diagonal matrix with the 
input data vector r[i), s*[i] is a DJM x 1 vector with the 
columns of To[i] stacked on top of each other and and the 
fl x 1 vector Vd contains a 1 in the rf-th position and zeros 
elsewhere. 

In order to analyze the proposed joint optimization proce- 
dure, we can rearrange the terms in x[i] and define a single 
D(JM + 1) x 1 parameter vector q[i] = [w T [i] s^[i]] T . We 
can therefore further express x[i] as 

x[i]=q H [i] ° Dxl °DxDJM 
[ 5R[i] Qd.jmx.djm 

= q H [i}U[Mi} 



dq H [i] dq\i) 

and testing if the terms are positive semi-definite. Specifically, 
H is positive semi-definite if a H Ha > for all nonzero a G 
c d(.jm+i)xD{jm+i) g2) 5 gyj Therefore, the optimization 

problem is convex if the Hessian H is positive semi-definite. 

Evaluating the partial differentiation in the expression given 
in d40b yields 

H = J2(Q H m[l}Q[^ d*[l])U[l] + £ U H [l]q[M^ H U[l} 

i=i i=i 

+ j^iq 11 mm d[i])u H [i] + £ u[i] q m^ H u H m 

i=i i=i 

(41) 

[*By examining H, we verify that the second and fourth terms 
are positive semi-definite, whereas the first and the third terms 
are indefinite. Therefore, the optimization problem can not be 
classified as convex. It is however important to remark that our 
studies indicate that there are no local minima and there exists 
multiple solutions (which are conjectured to be identical). 

In order to support this claim, we have checked the impact 
on the proposed algorithms of different initializations . This 
study confirmed that the algorithms are not subject to perfor- 
mance degradation due to the initialization although we have 
to bear in mind that the initialization Td{0) = Ojmxd anni- 
hilates the signal and must be avoided. We have also studied 
a particular case of the proposed scheme when JM = 1 and 
D = 1, which yields the cost function 



C(T D ,w) = E[\d-wT D r\ 2 }} 



(42) 



(38) 



By choosing Td (the "scalar" projection) fixed with D equal 
to 1, it is evident that the resulting function C(w, Td = 1, r) = 
\d — w* r| 2 ] is a convex one. In contrast to that, for a time- 
varying projection Td the plots of the function indicate that 
the function is no longer convex but it also does not exhibit 
local minima. The problem at hand can be generalized to the 
vector case, however, we can no longer verify the existence 
of local minima due to the multi-dimensional surface. This 
remains as an interesting open problem to be studied. 

B. Proof of Convergence of the Method 

In this subsection, we show that the proposed reduced- 
rank algorithm converges globally and exponentially to the 
optimal reduced-rank estimator Ifj"5l . ll2"2"l . ll2"3l . An issue that 
remains an open problem to be investigated is the transient 
behavior of the proposed method, which corresponds to the 
most significant difference between the proposed and existing 
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(MSWF and AVF) methods is on the transient performance. To Substituting d3T~b and using < f50b in d49l l, and manipulating the 
our knowledge, there exists no result for the transient analysis algebraic expressions, we can express (l49l in a more compact 
of the MSWF and the AVF methods, even though it has been way that is suitable for analysis, as given by 

w[i] = A 2 w[i-l}(w H [i-l]A 2 w[i-l)y 1 w H [i-l]w[i-i\. 

(52) 

(53) 



w[i] = Q[i] Q[i-1] ... Q[l] w[0], 
Q[i) = A 2l w[0}(w H [0}A 4l - 2 w[0})~ 1 w H [0}A 21 - 2 . 



(54) 



reported (and also verified in our studies) that the AVF 

has a superior convergence performance to the MSWF. 

As discussed in the previous subsection, the optimal solu- 

rr j - xx ,i i ■ i The above expression can be decomposed as follows 

tions To,opt an d w opt are not unique. However, the desired F F 

product of the optimal solutions, i.e., w opt = T n :0pt w ovt 

is known and given by RT 1/2 (R~ 1/2 p) VD 01,122), EI), 

where RT 1 ^ 2 is the square root of the input data co variance wnere 

matrix and the subscript 1 : D denotes truncation of the 

subspace. 

In order to proceed with our proof, let us rewrite the At this point, we need to establish that the norm of T D [i] for 
expressions in (JT3J and (O for time instant as follows all i is both lower and upper bounded, i.e., < ||Td[z]|| < 

oo for all i, and that w[i] = Tjj[i]w[i] approaches tu ptW 
-R[0]Td[0].R 1u [0] = f d[0] = p[0]w [0], (43) exponentially as i increases. Due to the linear mapping, the 

boundedness of To[i] is equivalent to that of w[i]. Therefore, 
we have upon convergence tu^ — 1] = w H [i — l]w[i — 1|. 
Since Hio^ijio^ — 1]|| < \\w[i — l]||||to[i]|| and ||io- H '[i — 
l]w[i — 1]|| = \\w[i 



(44) 



R[0]w[l] =T*[0}R[0}T D [0}w[l} = p[0], 
Using d43l we can obtain the following relation 

Rv>[0] = (TZ[0}R 2 [0}T D [0}y 1 T^ ) [0}R[0}p[0}w H [0l 

(45) 

Substituting the above result for R w [0] into the expression in 
d43l we get a recursive expression for Td [0] 

T D [0] = R[0]- 1 p[0] l I; H [0](Tg[0] J R[0]p[0]iI; i? [0])" 1 x 
x (t" [0}R 2 [0]T d [o]) -1 , 

Using (l44t we can express w[l] as 



l]w[i — 1] implies \\w[i 
\\ w [°°]\\ > \\ w [i 



2 , the relation w [i]w[i — 1] = 

> \\w\i — 111 | and hence 

> \\w 



(55) 



(46) 



In order to show that the upper bound ||td[oo]|| is finite, let 
us express the JAI x JM matrix Q[i] as a function of the 

JM x 1 vector w[i] = 
Ai 



Wi[l\ 

w 2 \i\ 



and the JM x JM matrix 



A = 



Substituting the previous expressions of 



w[i] and A into Q[i] given in d54b . we obtain 



w[l] = (TZ[0]R[0]T D [0]) l Tg[0]p[0], (47) 
Taking into account the relation w[l] = Td[1]w[1], we obtain 

W [i] = RiirM^miT^mmi^ii})- 1 - 



QH] 



A^w^O] 
A 2i w 2 [0] 



(«;f[0]Af- a «;i[0] 
W X[Q]At % - 2 wMr 



(56) 



(Tg [1}R 2 [1}T D [1]) \T% [0]R[0]T D [0]) 'Tf [0]p[0]Applying the matrix identity (A + B) 



(48) A~ 1 B)~ 1 A~ 1 to the decomposed Q[i] in 

More generally, we can express the proposed reduced-rank LS s ' s 
algorithm by the following recursion 



w[i] =T D [i]w[i] 

= R[i\- 1 p\i}w H [i](T^\i)R[i]p\i)w H \i})- 1 - 
■{T^R^Tnli])- 1 . 

■ (Tg[i - l}R[i - l]T D [i - l]) _1 Tg[i - l]p[i - 1]. 

(49) 

Since the optimal reduced-rank filter can be described by the 
SVD of R~ 1/2 p 02), G2,E1, where i?T 1/2 is the square 
root of the covariance matrix R and p is the cross-correlation 
vector, then we have 



Q[i] = diag(l_^ : l 0_^0,, 

D JM-D 



-A~ 1 B{I + 
and making i 

0(e[i]). (57) 



where e[i] = (A r+ i/A r ) 2i with A r +i and A r are the (r + l)th 
and the 7th largest singular values of R~ 1 ^ 2 p. From d57l i, it 
follows that for some positive constant k, we have < 
1 + ke[i]. From (l53l l. we obtain 



|iu[ooj|| < 
< 



(50) 



Considering that there exists some w[0] such that the randomly 
selected Td[0] can be written as 



(51) 



< 



T D [0] = R- 1/2 $w[0}. 



Q[oo]||...||Q[2]|| ||Q[1]|| \\Q[ 
11(1 + Mil) 



w 



w 



w 



OO 

exp(^; .g(l + /ce[z]) 



(58) 



exp 



cxp 



i=l 



(l 



k 



- {\ r+ i/\ r y 
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With the development above, the norm of w [i] is proven to be 
both lower and upper bounded. Once this is established, the 
expression in (|49b converges for large i to the reduced-rank 
Wiener filter. This can be verified by equating the terms of 
d52i i. which yields 



be obtained with an automatic rank adaptation algorithm for 
fine tuning the used rank. 



N=16, K=10 users, E IN =7 dB, f T=0.00005, Data support: 500 symbols 



W III 



- l p[i\w H 



l]R[i- l]T D [i- 
= R- 1/2 ® 1 A 1 <f>?p+ 0(e[i\). 



i]R[i\p[i\w 
1]) 



-•Hi 



]R 2 \i] 



i] P [i - 1] 



(59) 



where <J?i is a JM x D matrix with the D largest eigenvectors 
of R and Ai is a D x D matrix with the largest eigenvalues 
of R. 

VII. Simulations 

The performance of the proposed scheme and algorithms 
is assessed in terms of the uncoded bit error rate (BER) 
via simulations for space-time interference suppression in a 
DS-CDMA system. We consider dynamic fading situations, 
perfect synchronization and the proposed and existing adaptive 
algorithms are employed to adjust the filters and track the 
channel variations. Specifically, in our proposed reduced-rank 
STAP the output of the receiver x[i] is the input to a slicer 
that makes the decision about the transmitted symbol bk [i] for 
user k as follows 

h[i\ = Q{x\i}) = Q(w ff [i]Tf ]i]r\i]) (60) 

where Q ( • ) is the function that implements the slicer and the 
fcth user is assumed to be user 1. 

For all simulations, we use the initial values w[0) = 
[10 ... 0] T and T D [0] = [I D DJM -d} T - We assume 
L = 9 as an upper bound, employ QPSK symbols and 3- 
path channels with a power delay profile B31 given by 0, —3 
and —6 dB, where in each run the spacing between paths is 
obtained from a discrete uniform random variable between 1 
and 2 chips and the experiments are averaged over 200 runs. 
The power and the phase of each path is time-varying and 
follows Clarke's model |43l. This procedure corresponds to 
the generation of independent sequences of correlated unit 
power Rayleigh random variables for each path. The DoAs of 
the interferers and the desired user are uniformly distributed 
in (0, 2tt/3). The system has a power distribution among the 
users for each run that follows a log-normal distribution with 
associated standard deviation equal to 1.5 dB. We compare the 
proposed scheme with the Full-rank Oj], the MSWF |25] and 
the AVF 1 28 1 techniques for the design of linear space-time 
receivers and also the rank selection algorithms reported in 
ll25l and the ll32l with the proposed rank selection techniques. 

In the first scenario, we consider the BER performance 
versus the rank D with optimized parameters (forgetting 
factors A = 0.998) for all schemes. The results in Fig. [3] 
indicate that the best rank for the proposed scheme is D = 4 
for a data record of 500 symbols as it is very close to the 
optimal linear MMSE estimator. Studies with systems with 
different processing gains and loads show that D does not vary 
significantly with either the system size or the load. However, 
it should be remarked that considerable performance gains can 




Rank (D) 

Fig. 3. BER performance versus rank (D). 

In a second experiment, the BER convergence performance 
in a mobile communications situation is shown in Fig. 0] The 
channel coefficients are obtained with Clarke's model B31 
and the adaptive estimators of all methods are trained with 
200 symbols and are then switched to decision-directed mode. 
The results show that the proposed scheme has considerably 
better performance than the existing approaches and is able to 
adequately track the desired signal. In particular, the proposed 
reduced-rank algorithm converges in 100 symbols for the case 
of J = 1, in about 200 symbols for the case of J = 2 and 
in about 400 symbols for J = 3. This is substantially faster 
than the existing reduced-rank schemes, namely, the MSWF 
and the AVF (which are known to have the best performances 
available in the area) and the full-rank RLS algorithm. 



1 6, K=1 users, E./N =7 dB, f J=0.00005 




Number of received symbols 
Fig. 4. BER performance versus number of received symbols. 

In practice, the rank D can be adapted in order to obtain 
fast convergence and ensure good steady state performance 



10 



and tracking after convergence. To this end, we developed 
the automatic rank selection algorithm in Section V. We will 
assess this algorithm in a scenario identical to the previous 
experiment. The results in Fig. show that significant gains 
can be obtained from the use of the automatic rank selection al- 
gorithm. Specifically, we can notice that the proposed reduced- 
rank algorithm has a very fast convergence with D = 3 
even though it does not provide a steady state performance 
close to the full-rank optimal linear MMSE estimator. When 
the proposed reduced-rank algorithm employs D = 8 the 
convergence is notably slower even though it is able to 
approach the full-rank optimal linear MMSE estimator in 
steady state as shown in Fig. [5] and evidenced in our studies. 
Interestingly, when equipped with the proposed automatic rank 
selection algorithm the proposed reduced-rank RLS algorithm 
achieves a convergence performance as good as with D = 3 
and a steady state performance equivalent to that with D = 8. 
Another important issue is that the differences in performance 
are more pronounced for larger filters, when the usefulness of 
the automatic rank selection algorithm becomes more clear. 



N=16, K=10 users, E./N =7 dB, f T=0.00005 
bo d 




500 1000 1500 

Number of received symbols 

Fig. 6. BER performance versus number of received symbols with 
different automatic rank adaptation algorithms and the proposed 
reduced-rank scheme and algorithm. 




500 1000 1500 



Number of received symbols 

Fig. 5. BER performance versus number of received symbols with 
automatic rank adaptation. 

In order to assess the performance of the proposed rank 
selection algorithms, we consider the scenario of the previous 
experiment with J = 1,3 and compare the rank selection algo- 
rithms based on a stopping criterion ll25l . the cross-validation 
method of l32| and the proposed LS-based method with two 
variations, namely, the multiple filters and the extended filters 
approaches. The results shown in Fig. [6] indicate that the LS- 
based methods are slightly better than the other techniques. 
The cross-validation approach has the advantage that it does 
not require the setting of £> m i n and D max , however, it may 
perform a search over a higher range of values that leads 
to higher complexity. The remaining techniques operate with 
■Dmin = 3 and D max = 8. The method with a stopping rule 
has a performance slightly worse than the remaining schemes 
and its complexity is higher than the LS-based techniques due 
to the computation of the orthogonal projection. 

At this point, we will consider a study of the BER perfor- 



mance against the normalized fading rate of the channel (fdT) 
in the experiment shown in Fig. [7] We assess the performance 
of the receivers with a data record of 1000 symbols of training. 
The proposed algorithm is equipped with the automatic rank 
selection algorithm and the MSWF and the AVF algorithms are 
also equipped with the rank adaptation techniques reported in 
ll25l and Il32l . respectively. We observe from the curves in Fig. 
[7] that the proposed reduced-rank algorithm obtains substantial 
gains in BER performance over the existing MSWF and AVF 
algorithms and the full-rank RLS algorithm. We can notice 
that as the channel becomes more hostile the performance of 
the analyzed algorithms degrade, indicating that the adaptive 
techniques are encountering difficulties in dealing with the 
changing environment and interference. This behavior is more 
pronounced when the algorithms have to adjust filters with 
more coefficients, e.g. for more antenna elements (J = 2,3). 
In this regard, the reduced-rank algorithms obtain significant 
gains over the full-rank RLS algorithm and, in particular, 
the proposed reduced-rank algorithm achieves the best per- 
formance among them. 

The last experiment shows the BER performance versus 
the Eb/No and the number of users (K), which is illustrated 
in Fig. |8] In this scenario, all algorithms are trained with 
200 symbols and are switched to decision-directed mode for 
processing another 1500 symbols. The curves are obtained 
after 5000 runs. The proposed algorithm is equipped with the 
automatic rank selection algorithm and the MSWF and the 
AVF techniques are also equipped with the rank adaptation 
methods reported in ||251 and ll32l . respectively. The results 
confirm the excellent performance of the proposed reduced- 
rank algorithm, which can approach the performance of the 
optimal MMSE full-rank linear estimator (denoted simply as 
MMSE) that assumes the knowledge of the channels, the DoAs 
and the noise variance. In particular, the proposed reduced- 
rank algorithm can save up to 2 dB in E^/No in comparison 
with the existing reduced-rank techniques for the same BER 
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Fig. 7. BER performance versus number of received symbols. 



performance, whereas it can accommodate up to 4 more users 
than the MSWF and the AVF for the same BER performance. 
Interestingly, the performance of the optimal reduced-rank 
linear MMSE estimator [15] that assumes the knowledge of 
R and employs SVD is quite similar to the optimal full-rank 
one. For this reason, we only show the performance of the 
full-rank optimal linear MMSE estimator. 



(a) N=16, K=8 users, f T=0.0001 



(b) N=16, E IN =12 dB, f T=0.0001 




E h /N (dB) 



12 16 
Number of users (K) 



Fig. 8. BER performance against (a) E^/Nq (dB) and (b) Number of Users 
(K) for different techniques. 



VIII. Conclusions 

We proposed a reduced-rank scheme based on joint iterative 
optimization of parameter vectors. In the proposed scheme, 
the full-rank adaptive filters are responsible for estimating 
the subspace projection rather than the desired signal, which 
is estimated by a small reduced-rank filter. We developed a 
computationally efficient RLS algorithm for estimating the 
parameters of the proposed scheme and an automatic rank 
selection algorithm for computing the rank of the proposed 



RLS algorithm. The proposed algorithms do not require an 
SVD for dimensionality reduction and any knowledge about 
the order of the reduced-rank model. The results for space- 
time interference suppression in a DS-CDMA system show 
a performance significantly better than existing schemes and 
close to the full-rank optimal linear MMSE estimator in 
dynamic and hostile environments. The proposed algorithms 
can be employed in a variety of applications including spread 
spectrum and MIMO systems, wireless networks, cooperative 
communications and navigation receivers. 
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